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Discrete rearranging disordered patterns, part I: Robust 
statistical tools in two or three dimensions 
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Abstract. Discrete rearranging patterns include cellular patterns, for instance liquid foams, biological tis- 
sues, grains in polycrystals; assemblies of particles such as beads, granular materials, colloids, molecules, 
atoms; and interconnected networks. Such a pattern can be described as a list of links between neighbour- 
ing sites. Performing statistics on the links between neighbouring sites yields average quantities (hereafter 
"tools") as the result of direct measurements on images. These descriptive tools are flexible and suitable 
for various problems where quantitative measurements are required, whether in two or in three dimensions. 
Here, we present a coherent set of robust tools, in three steps. First, we revisit the definitions of three 
existing tools based on the texture matrix. Second, thanks to their more general definition, we embed these 
three tools in a self-consistent formalism, which includes three additional ones. Third, we show that the 
six tools together provide a direct correspondence between a small scale, where they quantify the discrete 
pattern's local distortion and rearrangements, and a large scale, where they help describe a material as a 
continuous medium. This enables to formulate elastic, plastic, fluid behaviours in a common, self-consistent 
modelling using continuous mechanics. Experiments, simulations and models can be expressed in the same 
language and directly compared. As an example, a companion paper [TJ provides an application to foam 
plasticity. 

PACS. 62.20.D Elasticity - 62.20.F Deformation and plasticity - 83.50.-v Deformation and flow in rheology 



1 Introduction 

Cellular patterns include: liquid foams or emulsions (Fig. 
[JJ); crystalline grains in polycrystals; or biological tissues 
(Fig. [2]). Assemblies of particles (Fig. [3]) include collec- 
tions of beads, molecules, or atoms; granular or colloidal 
materials; sets of tracers dispersed in a material, such as 
fluorescent probes or passively carried particles. Despite 
their tremendous diversity of sizes and physical properties, 
all these patterns have a common point: they are made 
of a large number of well-identified individual objects. 
We call them discrete patterns, where the word "discrete" 
here means the opposite of "continuous". Other discrete 
patterns include interconnected networks, of e.g. springs, 
polymers, biological macromolecules, fibers, or telecom- 
munication lines. 

We define the pattern as rearranging if the mutual ar- 
rangement of the individual objects can change. This is 
the case if they can move past each other, for instance due 
to mechanical strain (Figs. UJorEK), spontaneous motility 
(Fig. [lb), or thermal fluctuations (Fig. [3]). This is also the 
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case if the number of individual objects can change, for in- 
stance due to cell division or death (Fig. [lb), coalescence 
or nucleation of bubbles, shrinkage during coarsening of 
polycrystals or foams. 

Stimulated by the various imaging techniques, ref. [9] 
reviews many tools available to describe and quantita- 
tively characterise a pattern (that is, a single image). An- 
other tool is the texture: it appears in various contexts, in- 
cluding the order parameters of nematics, the microstruc- 
ture of polymers, or the fabric of grains; and has been used 
to describe mechanical strains by Aubouy et al. (see pjj] 
and references therein). It describes statistically how the 
individual objects are arranged with respect to each other. 
With a few simple measurements performed directly from 
an image, it extracts quantitative information relevant to 
the size and anisotropy of the pattern. Since it is based on 
statistics, it is particularly useful for disordered patterns. 
Here, we present a coherent set of robust tools (listed in 
Table [JJ, with a triple goal. 

First, we revisit existing definitions of the texture M, 
as well as the statistical strain U [10] and the rearrange- 
ments T [TT] based on it. In fact, cellular patterns are 
better characterised using cell centers than using their 
vertices. This remark enabled ref. [TT] to define a prelimi- 
nary version of T. Here, we show that it also yields a more 
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(b) i< 



Fig. 1. Liquid foams, (a) Heterogeneous flow: from left to 
right, around an obstacle [2]; liquid fraction ~ 10~ 4 , image 
width: 15 cm. (b) Homogeneous shear: in a rectangular box, 
deformed at constant area [3]; liquid fraction ~ 5 10 -2 , image 
width: 18 cm, courtesy C. Quilliet (Univ. Grenoble). 
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Table 1. Symmetric matrices used in the text. Equation num- 
bers correspond to their definitions. For comparison, the last 
row indicates the strains defined in continous mechanics for 
elastic, plastic and fluid behaviours. 



general definition of the texture, valid for all patterns. In 
addition, since a cell center is measured as an average over 
all pixels in a cell (while a vertex is a single pixel), mea- 
surements in experiments or simulations are more robust. 
The companion paper [TJ shows that displacements of cell 
centers (but not of cell vertices) are close to affine dis- 
placements. In order to make this paper self-contained, 




Fig. 2. Other cellular patterns, (a) Grains in a polycrystal of 
ice which rearranged during ice accumulation [4], image width: 
10 cm, courtesy J. Weiss (Univ. Grenoble), (b) Tissue of cells 
rearranging during the formation of a fruit fly (Drosophila) em- 
bryo fjj]: this thorax epithelium is labeled by the expression of 
the cell-cell adhesion molecule E-Cadherin-GFP; image width: 
160 /im, courtesy Y. Bellai'che (Inst. Curie). 



we recall and hopefully clarify the definitions of M and 
U. We also present a more general definition of T and 
derive explicitly its prefactor. 

Second, thanks to their more general definition, we em- 
bed these three tools in a self-consistent formalism, which 
includes three additional ones: B, V and P. From two 
successive images in a movie, we extract information re- 
garding the magnitude and direction of strain rate and 
rearrangements . 

Third, we show that the six tools together provide a 
direct correspondence between a small scale, where M, 
B and T quantify the discrete pattern's local distortion 
and rearrangements; and a large scale, where U, V and P 
help describe a material as a continuous medium without 
any details related with the discrete scale. This enables to 
formulate elastic, plastic, fluid behaviours in a common, 
self-consistent modelling using continuous mechanics even 
for a discrete material. Experiments, simulations and the- 
ories can be expressed in the same language to be directly 
compared. 

The only requirement is that the image should be of 
sufficient quality to extract the positions of the centers 
of each individual object (cell or particle); as well as the 



F. Graner et al.: Discrete rearranging patterns: Robust tools 



3 




Fig. 3. Assemblies of particles, (a) Beads repelling each other; 
they are placed on a vibrating loudspeaker, with an effect 
shown to be equivalent to thermal fluctuations [6]; image size: 
11.6 cm, courtesy G. Coupier (Univ. Grenoble) [7]. (b) Simu- 
lation of amorphous systems of atoms interacting via Lennard- 
Jones potential [8]: circles indicate each particle's effective ra- 
dius (here with a 20% dispersity), tangent circles correspond 
to vanishing interaction force; image size: arbitrary, courtesy 
A. Tanguy (Univ. Lyon). 



list of neighbour pairs (which objects are neighbours). 
All tools here are either static or kinematic, and rely on 
the image only; that is, they are independent of dynam- 
ics (stresses, masses and forces). They apply to discrete 
patterns regardless of the size of their individual objects, 
which can range from nanometers to meters or more. They 
regard simulations as well as experiments, and should en- 
able quantitative comparison between them. They apply 
whatever the pattern's disorder is. 

Our equations are valid in any dimensions. For clarity, 
we write them in 3D, and show that is is straightforward 
to rewrite them in 2D, see section [2.2.31 We specifically 
choose to illustrate this paper with 2D images (Figs.[T][3]), 
which are simpler and more common that 3D data. 

More precisely, we illustrate each definition on the ex- 
ample of a foam flow (Fig. [IJl) [2], which is both our orig- 
inal motivation and the most suitable example. Nitrogen 
is blown into water with commercial dishwashing liquid. 
Bubbles enter a channel, of length 1 m (only partly visi- 
ble on the picture), width 10 cm, and thickness 3.5 mm: 



a monolayer of bubbles forms (area Abubble = 16-0 mm 2 ), 
sandwiched between two glass plates (quasi-2D foam, liq- 
uid fraction less than a percent) . It steadily flows from left 
to right without vertical component (true 2D flow) until it 
reaches the free end of the channel. Coalescence and age- 
ing are below detection level. A 3 cm diameter obstacle is 
inserted into the foam channel. The foam is forced to flow 
around it, resulting in a spatially heterogeneous velocity 
field. Different regions simultaneously display different ve- 
locity gradients, internal strains, and rearrangement rates, 
and allow to sample simultaneously many different condi- 
tions. Bubbles naturally act as tracers of all relevant quan- 
tities; and on the other hand the foam's overall behaviour 
appears continuous. The total strain rate is partly used to 
deform bubbles and partly to make them move past each 
other; the companion paper [1] studies how it is shared 
between both contributions. 

Sections 130 and l2~2l start from the static description of 
Ref. [TO] and develop it step by step, for pedagogical pur- 
pose, while refining it. Section [231 describes the changes 
between two successive images. Section [3] is useful to com- 
pare measurements on different patterns; or to compare 
experiments and simulations. Section |4] is more theoreti- 
cal and regards specific applications: it discusses how to 
characterize materials which behave as continuous media, 
that is, where the quantities vary smoothly with space; 
and when it is possible to identify our statistical tools with 
the usual quantities of continuous mechanics. Appendices 
cover many practical, technical or theoretical details in- 
cluding all standard definitions and notations of matrices 
used in this paper. 

2 Texture and time evolution of links in the 
discrete pattern 

2.1 Ingredients 

The pattern is a collection of individual objects. Here we 
are interested in the relative positions of these objects, not 
in each object's shape (although both are related in partic- 
ular cases such as cellular patterns, see Appendix I A.3.3)) . 
We thus replace each object by a point called "site". 

The user should adapt the measurement tools to the 
pattern under consideration, and the scientific questions 
to be answered. For that purpose, the user should begin 
by deciding: (i) what are the relevant links, that is, pairs 
of sites which are connected (section 12 . 1 . 1|) ; and (ii) the 
averaging procedure (section 12.1.21 and Appendix IA.1[) . 

These choices are conventions, and thus rather free. 
The results of the measurements depends on the chosen 
definition, but they are much more robust than scalar 
measurements (see section 12". 3. 3p . Moreover, as long as the 
same definition is used for all measurements, the equations 
that relate the measurements of the different quantities 
(such as eq. [9]) are valid independently of the chosen defi- 
nition. Once conventions are chosen, it is thus important 
to use them consistently. 
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Fig. 4. Definition of sites and links, (a) Cellular pattern. Back- 
ground: detail from Fig. Jib)- Foreground: a site is a cell's geo- 
metrical center; there is a link between two centers if their cells 
touch, (b) Particle assembly. Background: detail from Fig. J3K). 
Foreground: a site is a particle's position; the links are defined 
as discussed in the text (here a Delaunay triangulation). 

2.1.1 Links between neighbouring sites 

In a cellular pattern (Figs. Ill2p . it is often advisable to 
choose as sites each cell's geometrical center (see Fig. [4^). 
However, alternative choices exist. For instance, a user 
interested in studies of dynamics might prefer the center of 
mass, if different from the geometrical center. Similarly, a 
biologist might be more interested in the cell's centrosome 
or nucleus. Note that we do not advise to use a definition 
based on vertices (see Appendix IA.3,2j) , 

When two cells touch each other it defines that their 
sites are connected. This is unambiguous if cells walls are 
thin. This is the case for grains in polycristals, cells in an 
epithelium, or in a foam with low amount of water (Figs. 

If cells walls are thick, as is the case in a foam with a 
higher amount of water (Fig. Q})), different definitions of 
neighbours are possible. For instance, two cells are defined 
as neighbours if their distance is smaller than a given cut- 
off. Or, if they are neighbours on a skeletonized image; 
that is, after an image analysis software has reduced cell 
walls to one pixel thick black lines on a white background. 
If cell walls are too thick, cells are really separated (as in a 
bubbly liquid, where bubbles are round and far from each 
other) and can be treated like the particles, which we now 
discuss. 

If each object is a particle (as in Fig. [3]) it is natural to 
choose its center as site (see Fig. [4ja) . There are various 
possible choices for the links. Since the tools characterize 
patterns and not forces, the definition of links is indepen- 
dent of interactions: a link between sites does not mean 
that sites interact; conversely, sites which interact are not 
necessarily linked. Whatever the chosen definition, it is 
important that each particle has only a finite number of 
neighbours. 

In a first case (Fig. [3}d), the average distance between 
particles is comparable to their average radius; for in- 
stance, for a dense (also called compact or jammed) colloid 
or granular material. We then recommend to define that 
two particles are linked if their distance is less than a cho- 



sen cut-off. For hard spheres, this cut-off should be the 
sphere's radius plus a small tolerance. 

In the opposite case, the average distance between par- 
ticles is much larger than their average radius (Fig. [3^) ; 
for instance, for a decompacted colloid or granular ma- 
terial. We then recommend to recreate a cellular pattern 
by attributing to each particle its Voronoi domain (the 
set of points surrounding this particle, closer to it than to 
any other particle). One then chooses to define that two 
particles are linked if their Voronoi domains touch; this is 
called the "Delaunay triangulation" of the particles. 

If the pattern is a network, it is natural to choose the 
nodes as sites. The connexions are physically materialised, 
and thus unambiguously defined. 



2.1.2 Averaging 

The present tools aim at describing the collective proper- 
ties of links. In what follows, (.) denotes the average over 
a set of links relevant to the user: 

(•> = ^£(-). 

where the sum is taken over the number Ntot of such links. 
Appendix I A. II presents some technical details, especially 
regarding the boundaries of the averaging region, which 
can be treated as sharp or smooth. 

The scale of study determines the number of links in- 
cluded. Performing the same analysis at different scales 
(Fig. [6j) enables to obtain multi-scale results |4 |fTT1fT2) . For 
instance, we can measure the dependence of pattern fluc- 
tuations with scale for particle assemblies. 

Choosing to average over a small number Ntot of links 
yields access to detailed local information. For instance, 
the local heterogeneity of a sample of ice can be measured 
by including the links around one single grain, then per- 
forming a comparison between different grains [4j. Sim- 
ilarly, to study the anisotropy of a biological cell which 
divides, one can include only the links starting at this 
cell's center [5]. 

On the other hand, choosing a large number Ntot of 
links yields better statistics. This is the case for instance 
if the system is homogenous in space. In a homogenously 
sheared foam (Fig. QJd), it makes sense to consider that 
all bubbles play a similar role, and average over the whole 
foam. Averaging over all links contained in the whole im- 
age enables to detect the overall anisotropy of an ice sam- 
ple or an epithelium, and compare it with other samples 

BE] • 

Even if the system is invariant along only one direction 
of space, one can average over this direction [12]. Similarly, 
in a flow which is invariant in time, one can average over 
time (Fig.QJ,). In what follows, figures are prepared with 
1000 successive video images, representing millions of bub- 
bles: a time average yields good statistics and details of 
local variations, even if only a small part (i.e. few links) 
of each image is included (Fig. [6fc) . 
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2.2 Texture M: current state of the pattern. 

We include this section, already published [lOlfTT] . in order 
to make the present paper self-contained. 



2.2.1 Definition and measurement 

A pair of neighbour sites of coordinates rx = (xx,yx, Zi) 
and T2 = {x2, y% , Z2) constitutes a link. Reducing the pat- 
tern to a set of links sets aside the detailed information 
regarding the actual positions of each site. 
The link vector: 



r 2 - rx, 



(1) 



has coordinates (X, Y, Z) = (x-i — Xx, J/2 — Hi, %i — #1). It 
carries the information on link length and angle. However, 
i and — £ play the same physical role: an average over sev- 
eral ^s will yield a result which depends on this arbitrary 
choice of sign (and, in practice, if there are enough links, 
the average (I) turns out to be close to zero). 

On the other hand, the number £ 2 = X 2 +Y 2 + Z 2 is in- 
variant under the change £—* —£ and thus has a physically 
relevant (and non-zero) average: (^ 2 ) = (A 2 + Y 2 + Z 2 }. 
It reflects the average square link length, but loses the 
information of angle. 

The link matrix m combines the advantages of both: 



X 2 XY XZ\ 
YX Y 2 YZ 
ZX ZY Z 2 J 



(2) 



Its trace is Tr (m) 

m 

M = (rri 



Its average defines the texture 



(X 2 ) (XY) (XZ)\ 
(YX) (Y 2 ) (YZ) 
(ZX) (ZY) (Z 2 ) J 



(3) 



It is expressed in m 2 . As required, it stores the same infor- 
mation regarding the current pattern: the square length, 
readily visible as the sum of diagonal terms; the angle and 
magnitude of anisotropy, as discussed below. 



2.2.2 Diagonalisation and representation 





Fig. 5. Measurement of texture. Snapshots of two regions se- 
lected in Fig. (JT^i): the foam is nearly isotropic in A, not in 
B. From the statistical analysis of links (lines), and time aver- 
age over several images, we compute the corresponding M . We 
represent it by an ellipse with axes proportional to the eigen- 
values \i\ in A it is nearly circular. Thin lines indicate the axes 
with positive eigenvalues (i.e. here all axes). 



M can be represented as an ellipsoid, which axes di- 
rections are that in which M is diagonal, represented as 
thin lines on Fig. |(5]). Each ellipsoid's axes length is pro- 
portional to the corresponding A^. It is expressed in m 2 : 
this is less intuitive than m, and ellipsoids are more elon- 
gated that the actual cell shape; but this is necessary for 
the consistency with the representation of the other ma- 
trices (see Appendix lA. 3.3ft . The square link length (f 2 ) is 
reflected in the size of the ellipsoid, more precisely as the 
square root of the sum of the three axes lengths; it is thus 
not proportional to the ellipsoid's volume. The direction 
in which links are longer is represented by the direction of 
ellipsoid elongation: the greater the pattern's anisotropy, 
the more elongated the ellipsoid. If the texture is measured 
at several regions of the image, it is represented as several 
ellipsoids, that is, a map of the texture field M(R, t) (see 
also section [2331 and Fig. [6]). 

When the pattern is statistically isotropic, so is its 
texture. It is thus diagonal in any system of axes, and the 
three AjS are exactly equal, A^ = (£ 2 ) /3: 



M 



isotropic 




(5) 



By construction, M is a matrix with symmetric off-diagonal That is, the texture of an isotropic pattern contains only 
terms (XY = YX, etc.). It can thus be diagonalised (see 
Appendix IB. 2l for details): 



the information of length: M = (€ 2 ) J 3 /3, where I3 is 
the identity matrix in 3D. It is represented as a sphere. In 
that case, all axes are equivalent (or "degenerated"). 



diag M = 




(4) 



Its three eigenvalues A; (i = 1, 2 or 3) are positive. Their 
sum, TrAf, is exactly (^ 2 ). In practice, they usually have 
the same order of magnitude, of order of (^ 2 ) /3. 

In a truly 3D pattern, M has strictly positive eigen- 
values (except in unphysical examples). Thus its inverse 
M^ 1 always exists (eq. 



2.2.3 Two-dimensional case 

If the pattern under consideration is contained in a plane, 
as are most experimental images, we turn to a 2D nota- 
tion. As mentioned, this is straightforward: 



2D : 



_ / (A 2 ) (XY) \ 
M -{(YX) (Y 2 ))- 



(C) 
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There exist two orthogonal axes (eigenvectors) in which 
M would be diagonal: 



2D : diag M = 



Ai 
A 2 



(7) 



with strictly positive A„ (i = 1 or 2). If we call (£ + ) the 
r.m.s. length of links in the direction of elongation (say, 
1) and (£-) the r.m.s. length of links in the direction of 
compression (say, 2), then Ai ~ /2 and A2 ~ (jPL) /2. 
Its inverse M~ l always exists. 

In 2D, M is represented by an ellipse (Fig. HJ. Mea- 
surements can be performed at larger scale to decrease the 
noise due to fluctuations, or at smaller scale to evidence 
more details of the spatial variations (Fig. [6]) . 

When the pattern is isotropic, its texture is diagonal 
with any choice of axes: 



2D 
isotropic 







2 

v 1 



(8) 



That is, M = (i? 2 ) J 2 /2, where 1 2 is the identity matrix 
in 2D. All axes are equivalent (or "degenerated") and the 
angle of eigenvectors is not defined. M is represented by a 
circle and the thin lines lose their signification (Fig. UK). 



2.3 Time evolution 

Differentiating eq. determines how M varies. Appendix 
IC.1[ useful for practical calculations, discusses finite size 
effects due to the time interval between successive images 
of a movie. Neglecting these effects in eq. ([66)1 yields the 
simplified time evolution of M: 



DM 



M 



d log N t , 

at 



= - V -Jm + B 



(9) 



The variation of N to t is negligible in most physical ex- 
amples; however, it is significant for instance in biological 
tissues with many divisions [5] or in coarsening systems 
such as ageing foams p~3|TT4] . 

The three terms of the r.h.s. can be measured on a 
movie, and have the following meaning. In the time inter- 
val between two successive images, some links enter or exit 
the region of averaging; some links change their length or 
angle; some links are created or destroyed, respectively. 
They are now discussed one by one (sections 12.3.1) 12.3.21 
and 12.3.31 respectively). 
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Fig. 6. Map of the texture measured in each region of the foam 
flow (Fig. [Ik). The area Vbox of each averaging box corresponds 
to (a) 3, (b) 1 and (c) 0.3 bubbles (sub-bubble resolution); 
since there are 1000 movie images, and 3 times more links 
than bubbles [13], this corresponds to averages over 10 4 , 3 10 3 
and 10 3 links, respectively. Scale: actual image size, 15 cm x 
10 cm; M (ellipses): bar = 10 mm 2 , for the eigenvalues, which 
are here all positive, and represented as ellipse axes lengths. 



2.3.1 Flux J M 

In eq. ((9}, Jm is the flux of advection, that is, the trans- 
port of texture. It counts the rate at which links enter or 
exit throught the sides of the region of averaging. Tech- 
nically, it is a rank-three tensor (i.e. with 3 indices): for 



more details see ref. [15] and Appendix IB. 31 In a good 
approximation, Jm — v ® M, where v is the local av- 
erage velocity, see Appendix IC.lt Its divergence V • Jm 
counts the net balance between links that enter and exit; 
it vanishes if M is spatially homogeneous, or at least is 
constant along a flux line; it also vanishes if the local av- 
erage velocity is zero. 
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2.3.2 Geometrical texture changes: B 




Fig. 7. Changes in the shape of links (geometry). In a first 
snapshot of a small region of the foam, links are represented by 
dashed lines. The same bubbles are tracked on a next image, 
with links represented by solid lines. To evidence how each link 
changes, we have removed the overall translation (which plays 
no role here), and superimposed both snapshots. We calculate 
B and plot it as a coffee bean: an ellipse with a thin line 
indicating the positive eigenvalue (direction of extension). 

B describes the changes in the pattern's overall shape, 
that is, geometry: at which rate, and in which direction, 
the pattern deforms. It reflects relative movements: it is 
insensitive to a global, collective translation. 

B is based on the links which, on both successive im- 
ages, exist {i.e. do not undergo topological rerrangement) 
and belong to the region of averaging (i.e. are not ad- 
verted). These links may change in length and direction 
(Fig.[7J). We obtain each link's contribution directly from 
eq. ([2]) and average it, like in eq. Q: 

„ / dm \ , . 

B = \-dr)- (10) 

Appendix IC.ll provides more details; it also defines C, 
which symmetrical part is related with B, and will turn 
useful to define W and V in section l3.2.1l The right hand 
side of eq. (fit)]) should not be confused with d(m) /dt; 
that is, B is not equal to the variation of M: as eq. Q 
shows, the difference between them is due to changes in 
the links included in the averages. 

B is symmetric. Its units are in m 2 s _1 . It has a positive 
eigenvalue in a direction of extension, and a negative one 
in a direction of compression. In case of dilation, all its 
eigenvalues are positive. In a region of shear, it is plotted 
as an ellipse with one thin line drawn on it, similar to 
a "coffee bean" (Fig. [7|. We can plot a map of B: it is 
similar to Fig. i[T0|). data not shown. 



2.3.3 Topological texture changes: T 

T reflects the topological changes, namely, changes in the 
list of links: creation and destruction, that is, source term 
of the texture. Appendix IC.ll provides details. Briefly, each 
link t a which has appeared since the preceding image has 




Fig. 8. Changes in the list of links (topology). They are here 
illustrated by cellular patterns but apply to all other patterns 
as well, (a) Neighbour exchange in 2D: snapshots extracted 
from a dry foam pQ, one link disappears (dashes) and another 
appears (thick grey line); and corresponding representation of 
T as an ellipse, where a thin line indicates the positive eigen- 
value, (b) Same, sketched in 3D: three links (hence 3 faces) 
disappear and one appears, (c) Site disparition, in 2D or 3D: 
all its links disappear, (d) Coalescence of two sites, in 2D or 3D: 
the link between them disappears, the links to their common 
neighbours merge (here there is a total of five disappearances 
and two creations). 
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a contribution given by eq. ([2]), noted m a ; and similarly 
the contribution of a link which disappears before the next 
image is noted mj, Averaging over all links which appear 
or disappear between successive images defines T as: 



T = h a (m a ) - h d (m d ) . (11) 

The quantity n a (resp. rid), expressed in s -1 , is not the 
time derivative of a physical quantity (which we would 
note d/dt). It is the rate of link appearance (resp. dis- 
appearance), per unit time and per existing link. If n a 
and rid are equal, their inverse is the average link's life 
expectancy. 

T is expressed in m 2 s _1 . It characterises the total ef- 
fect on the pattern of all topological changes occurring 
between two images. By construction it is symmetric, like 
M: it can thus be diagonalized and represented as an ellip- 
soid. It is robust to artefacts and errors in determinations 
of neighbours |llj . It is general, and includes information 
of frequency, size, direction and anisotropy for all con- 
tributions of all topological changes: they can be treated 
indifferently and added together. However, as we now dis- 
cuss, the user might be interested in studying separately 
the contributions of the different processes [5]. 

The coalescence of two sites (Fig. [8ji) corresponds in 
foams to the breakage of a liquid wall between two bub- 
bles, with a net balance of minus one site [16]. The reverse 
process corresponds in epithelia to a cell division, and re- 
sults in one more site [17]. When the number of sites de- 
creases (resp. increases), so does the number of links, and 
T usually has only negative (resp. positive) eigenvalues. 
The variation in the number of sites and links is thus vis- 
ible in the trace of T. 

Neighbour exchanges (Fig. HK,b), also called "Tl" in 
the case of cellular patterns [HUH], preserve the sites; in 
2D, they also preserve the number of links, n a = fid- It 
is for that case that ref. [TT] introduced a specific defini- 
tion of T. T usually is mostly deviatoric (Appendix IB.ip . 
with both positive and negative eigenvalues (and a van- 
ishing trace), in directions correlated with the appearing 
and disappearing links. Note that the eigenvectors are ex- 
actly orthogonal, while the appearing and disappearing 
links need not be: thus eigenvectors are not strictly paral- 
lel to links, especially when T is measured as an average 
over several individual topological processes. For our flow- 
ing foam example, a map of T (data not shown) would be 
similar to Fig. 18 of ref. [TT], or to Fig. [TT] below. 

If a pattern has a free surface, when sites exchange 
neighbours the number of links can vary, and thus locally 
T can be far from deviatoric. 

The disparition of a site (Fig. [St) corresponds in foams 
to a bubble which shrinks, also called "T2" [T3lll4| : and in 
epithelia to a cell which dies, or exits the epithelium plane. 
The reverse process is a site nucleation. Both processes 
have an approximately isotropic contribution to T. 



3 Statistical tools to obtain relative 
deformations and their time evolution 

This section facilitates comparison between different ex- 
periments; or between experiments, simulations and the- 
ory. This is useful for large scale deformations and flows: e. 
g. of particle assemblies, of foams and emulsions, of gran- 
ular materials, or of biological tissues during morphogen- 
esis. 

Here we try to link the discrete, local description of 
Section [2] with the continuous, global description of Sec- 
tional For this continuous description to be self-consistent, 
it is necessary to get rid of the discrete objects' length 
scale, that is, the typical size of links. For each of the three 
discrete quantities M, £?, T defined in section[2j it is pos- 
sible to construct a continuous counterpart: that is, a tool 
which is dimensionless (or expressed in s _1 ), with no m 2 
any longer. This defines U, V, P, respectively (sections 
[3 HI3XT1 and 13X21 respectively). 



3.1 Statistical internal strain : U 

The internal strain has been defined by Aubouy et al. 
[TO] through a comparison between the current pattern 
and a reference one. We include it here in order to make 
the present paper self-contained, and to provide additional 
explanations and examples. 



3.1.1 Strain of a single link 

Consider first a link £ of length £, and apply to it an in- 
finitesimal variation dt. Its relative extension, or infinites- 
imal strain, is d£/£, or equivalently d(log^) [18]. The "true 
strain" (also called "Hencky strain" [18]) is defined with 
respect to a state £q chosen as a reference (often a state 
without stress) using several equivalent expressions: 

= i[logTr(m)-logTr(m )]. (12) 

We perform these manipulations because the last expres- 
sion of eq. (fT2"l) is the easiest to generalise. It is not a 
problem to take the log of dimensioned quantities (here, 
the square of a length) because this cancels out in the final 
result. 



3.1.2 Statistical strain of the pattern 

For a whole pattern, replacing Tr (m.) by M enables to 
perform statistical averages over links. The logarithm of 
M is unambiguously defined and is easily performed in 
three standard steps on a computer (Appendix IB.2|) . It 
suffices to first, switch to the three orthogonal axes (M's 
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eigenvectors) in which M is diagonal; second, take the log- 
arithm of its eigenvalues, which are strictly positive (sec- 
tion 



diag log (M) 




(13) 



and third, switch back to the original axes. It is neces- 
sary to perform first all linear operations such as aver- 
aging. This ensures in particular that all A^s in eq. I|13p 
are non-zero. Taking the logarithm, which is a non-linear 
operation, has to be performed later. 

Eq. fl3|) . like eq. fl"2"|) . requires to define a reference, 
expressed in the same units as M, so that the difference 
of their logarithms is well defined and dimensionless. Such 
a reference texture Mo is discussed in section 13.1.31 

The "statistical internal strain" is defined [101 as: 



U 



\{\ogM 



log M ) 



(14) 



Here U completely characterises the material's current 
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Fig. 9. Elastic behaviour. Map of the statistical internal strain 
U f eq. n~4|) measured on the foam of Fig. QJi) . Coffee bean axes 
indicate the direction and amplitude of stretching (indicated 
by a thin line) and compression. M is chosen as the averaged 
value of M measured at the left and right of the image, far 
from the obstacle. Same box size as in Fig. J6}j). Scale: for 
ellipses axes lengths, bar = 1 (dimensionless) for the positive 
eigenvalue and the absolute value of the negative one. 



strain: relative dilation, amplitude and direction of anisotropy. 



3.1.3 Reference texture Mo 



then M = fc 2 M . Thus U = \og(k)Io, as is expected for 
instance for gases; that is, in 3D: 



Practical details regarding the reference texture Mo are 
presented in appendix I A.21 

Eq. (fl4|) shows that the exact choice of M affects 
the value of U but not its variations. It thus does not 
appear explicitly in the kinematics feqs. I21|80|) nor in the 
dynamics (for instance in the value of the shear modulus, 
appendix lA.3.3j) . Moreover, eq. (fl4|) remains unchanged if 
we multiply both Mo and M by a prefactor; this is why 
the exact unit (e.g. m 2 , mm 2 , /im 2 ) in which Mo and M 
are expressed is unimportant, as long as it is the same unit 
for both. 

Whatever the choice, the reference is defined by the 
texture Mo- It suffices to determine 6 numbers (3 num- 
bers, if in 2D); or 1, in the (most common) case where 
M is isotropic. Only the reference texture corresponding 
to the current state plays a role; past changes of the ref- 
erence pattern, for instance during an irreversible strain 
(also called "work hardening" [IS]), need not be taken into 
account. It is never necessary to know the details of the 
corresponding pattern's structure, such as the positions of 
each object one by one: it is even not necessary that this 
pattern exists and is realisable. 



3.1.4 Examples 

This section presents a few examples and particular cases 
of internal strain. 

(i) If the material is uniformly dilated (affine defor- 
mation, see section 14.3. 1|) by a factor k in all directions, 



log/c 
logfc 
logfc 



(ii) Conversely, if the material is uniformly dilated by 
a factor fc in one direction and compressed by a factor 1/fc 
in another direction, then: 



diag U 



logfc 0' 
- log fc 




(iii) For incompressible materials, C/'s diagonal terms 
are usually both positive and negative, and their sum 
is usually small: U is mostly deviatoric (Appendix IB.ip . 
Note that even in incompressible materials the links' mean 
square length can vary slightly, so that Tr(C7) is not nec- 
essarily strictly zero. For instance, it reaches 0.03 in Fig. 
(JTJd) where a foam is sheared while keeping bubble num- 
ber and total foam area exactly constant (Ataei Talebi 
and Quilliet, private communication). 

(iv) In Fig. ([9}, most ellipses look circular; deviations 
from circles occur close to the obstacle. We distinguish re- 
gions where extension dominates, and ellipses are stretched 
like coffee beans, from regions where compression domi- 
nates, where the ellipses are flattened like capsules. 

(v) If (but only if) Mo is isotropic, then U is diagonal 
in the same axes as M. Then eq. (fl3|) enables to rewrite 
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eq. (fl4|) more explicitly: 



diag U 



\ 











(15) 



where Ao is Af 's eigenvalue (e.g. Ao = (£q) /3 is M if 
we use the definition of eq. |44|). Eq. I|15p reflects that M 
and U have the same eigenvectors: they commute. Eq. 
lfl5|) also relates the trace of U with Af's determinant 
(product of eigenvalues) : 



Tr U = log , 



' A1A2A3 



4 



- log (dot M) — - log (det M ) 



(16) 



(vi) in the limit of small strains, i.e. when M remains 
close enough to M , eq. lfl4|) can be linearised [10]. The 
difference of logarithms simply amounts to a division by 
M , that is: U ~ (M-M )M^ 1 /2. This is true whether 
M is isotropic or not (unlike eq. [15]) . This approximation 
is used in Appendix IC.2. 21 



3.2 Kinematics: time evolution 

3.2.1 Statistical velocity gradient: W and V 

We want to define the continuous counterpart of the geo- 
metrical changes B (eq.[T0|). We use Af _1 (eq. 155)1 . which 
is in m -2 , and is always defined. For reasons which appear 
below feqs. I271I30]) . we use C (eq. [71]) as an intermediate 
step, and define W as: 



W = M 1 C 



dt\ 

~dt) 



(17) 



W has the dimension of a strain rate (s^ 1 ): its order of 
magnitude is the links' average variation rate. Like B, it 
vanishes when the pattern moves as a whole, with a rigid 
body translation. 

For reasons which appear below (eqs. I271I30]) . we call 
it the "statistical velocity gradient", and we purposedly 
define it as iW _1 C rather than C M.~ x . In general, G is 
not symmetric. 

In practice, the most useful quantity is its symmetric 
part, the "statistical symmetrised velocity gradient": 



V = 



W + W* M~ l C + C t M~ 



(18) 



It is the rate of variation of U due to the links' stretch- 
ing and relaxation. In cases where C and M commute, 
M~ 1 C is symmetric and eq. lfT8]) simply writes V = 
BM~ 1 /2. 

Fig. (flQj) plots an example of V. It is large all around 
the obstacle, but only very close to it; it is almost the 
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Fig. 10. Fluid behaviour. Map of the statistical symmetrised 
velocity gradient V (eq. 1 1 8 measured on the foam of Fig. 
(Jlp,). Coffee bean axes indicate the direction and amplitude of 
stretching rate (indicated by a thin line) and compression rate. 
Same box size as in Fig. Et). Scale: for ellipses axes lengths, 
bar = 0.1 s _1 for the positive eigenvalue and the absolute value 
of the negative one. Grey levels: statistical vorticity from the 
rotation rate fl (eq. IT9jl in s _1 . 



same before and after the obstacle. When the material's 
density is constant, TrV is small (but not necessarily ex- 
actly zero), and the corresponding ellipse is nearly (but 
not necessarily exactly) circular. 

The anti-symmetric part is the statistical rotation rate: 



fl 



W -W l M- X C -C l M 



(19) 



It has 3 independent components in 3D, but only 1 in 2D 
(appendix IB. 1|) . Thus Fig. (fTD]l plots it as grey levels. 



3.2.2 Statistical topological rearrangement rate: P 

We define the continuous counterpart of the topological 
changes T (eq. HPT) in a way similar to eq. fT8]) : 



1 M^T + TM" 

2 2 



(20) 



Here we have introduced a factor —1/2 so that P is the 
term which unloads the statistical internal strain, as will 
appear in eqs. (|21I80[) . In cases where T and M commute, 
such as in the companion paper [Tj, M~ X T is symmetric 
and eq. J2Q]) simply writes P = —TM~ 1 /2. 

This "statistical topological rearrangement rate" P 
(eq. [2Q]) has the dimension of s _1 . It measures the fre- 
quency and direction of rearrangements: it is of the order 
of magnitude of the number of changes per unit time and 
per link. Corresponding ellipses are elongated like coffee 
beans (resp: flattened like capsules) if the number of links 
decreases (resp: increases); if the number of links is con- 
served, ellipses are nearly circular. 

As an example, Fig. (fTT]l shows that the rearrange- 
ments are more frequent just in front of the obstacle, or 
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Fig. 11. Plastic behaviour. Map of the topological strain rate 
P feq. l20jl measured on the foam of Fig. ([TJi). Coffee bean axes 
indicate the direction of links which have just disappeared (in- 
dicated by a thin line) and just appeared; note that this is the 
inverse of Fig. (|8|), due to the minus sign in eq. (| 20 j) . Same box 
size as in Fig. [6K). Measurement boxes touching the obstacle 
were removed. Scale: for ellipses axes lengths, bar = 0.1 s _1 for 
the positive eigenvalue and the absolute value of the negative 
one, both proportional to the frequency of rearrangements. 

in a very narrow region behind it. The rate of rearrange- 
ments decreases smoothly with the distance to the ob- 
stacle. This is due to the foam's elasticity. It contrasts 
with the sharp transition between solid-like and fluid-like 
regions observed in purely visco-plastic materials |20p21j. 
The companion paper pQ presents an example with a larger 
spatial distribution of topological events, which enables for 
a better spatial resolution. 



3.2.3 Kinematic equation of evolution 

We have thus three independent symmetric matrices: U, 
V and P. As discussed in Appendix IC.2.3[ there is a re- 
lation between them. 

In the case where we can neglect the variation in N tot 
and the higher order terms in U, the time evolution of M 
approximately simplifies as (eqs. I9I80|81[) : 

v = ^ + r. 

That is, the (statistical) symmetrised velocity gradient is 
shared between two contributions: one part (which in- 
cludes advection and rotation) changes the (statistical) 
internal strain, the other part is the (statistical) topolog- 
ical rearrangement rate. How V is shared between both 
contributions constitutes the main subject of the compan- 
ion paper [I]. 

Physically, eq. (j2Tj) means that, when a perturbation 
is applied to the overall shape of the pattern, part of it 
affects the appearance of the pattern (loading) and the 
other part goes into rearrangements (unloading). Section 
2] introduces a parallel point of view, in particular with eq. 



4 Continuous mechanics 

This section applies to materials which dynamics can be 
described using continuous mechanics (section PTTj) . in terms 
of stresses. We examine whether it is possible to relate 
the continuous, large-scale, dynamical description on one 
hand; and on the other hand the statistical measurements 
based on discrete objects introduced in section O which 
describe the pattern's connections (topology), shape (ge- 
ometry) and movements (kinematics). 

We recall that continous mechanics involves three kine- 
matical quantities e to t, £ei and t p i (section |4~2"1) . which are 
related through eq. 11251) . We then try (section l4~3j) to iden- 
tify it with eq. l[2"T]). 



4.1 Continuous description and RVE 

If the material acts as a continuous medium |19 [ [22 j [23]. 
it usually has the following properties. First, there exists 
a range of Vbox sizes over which measurements yield the 
same results In that case, the box is called a repre- 
sentative volume element (RVE) . This is usually obeyed if 
Vbox is much larger than the range of interaction between 
individual objects, and also larger than the correlation 
length of their disorder (but these conditions are neither 
necessary nor sufficient). Second, its description can be 
local in space, that is, its equation of evolution involves 
partial space derivatives, and the spatial variations of its 
solutions look smooth. Third, the average quantities have 
at large scale a role more important than that of fluctua- 
tions. 

Regarding the choice of the RVE, the discussion of 
section [2.1.21 applies. Here again, averages (.) on detailed 
geometrical quantities are performed on a spatial box of 
volume Vbox and over a time r selected to suit the problem 
under consideration. The shape of the box should prefer- 
ably respect the system's symmetries. 

For the present purpose of a continuous description, 
there is however the additional requirement that Ntot ^ 1- 

— 1/2 

More precisely, the relative statistical uncertainty N tot ' 
should be smaller than the relative precision required by 
the user. A few tens or hundreds of links are often enough 
(there is no need for 10 23 links). This does not set any 
theoretical lower limit to the size of Vbox '■ it can well be as 
small as the link size, or even smaller, if there are enough 
images to average (Fig. [6b). 



4.2 Elastic, plastic, fluid behaviours 

If the pattern behaves as a continuous material, we can 
consider a RVE (section |4~T]) at position R. The velocity 
field is (v) (R), that is, an average over the whole RVE. 
If i?i and i?2 are the positions of two RVEs, the velocity 
gradient Vd is the spatial derivative of the velocity field, 
and Vi>' is its transposed (cqs. I60|6"1~T) . then: 

(v) (R 2 ) =i (v) (ili) + Vt/ ■ (R 2 - Ri) . (22) 
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Details on this notation can be found in Appendix | 
Eq. (f22|) neglects terms of order of |i?2 — Ri | 2 and higher. 
It describes the velocity field as continuous and affine, 
that is, a term which varies linearly with position plus a 
constant term (offset). 

One of the key ingredients of continuous mechanics is 
the velocity gradient's symmetrical part, that is, the total 
strain rate: 



Vd + Vv' 



(23) 



This is a purely kinematical quantity, but it determines 
the contribution to the viscous (dissipative) stress [22] . 

For small strain (linear elastic regime), neglecting ad- 
vection and rotation, the integration of eq. ([23|l defines a 
total applied strain, which is a function of the past history 
of the sample, as: 



Etot 



dt Et 



Vu + Vu 



(24) 



Here Vu is the gradient of the displacement field u, and 
Etot its symmetrical part. 

The total strain rate Etot contributes in part (loading) 
to change the elastic strain e e ; , and in part (unloading) to 
a plastic strain rate £ p i which is defined by their difference: 



Stot 



T>£el 

vt 



' pi • 



(25) 



Fig. 12. (a) Affine and (b) non affine flow field. 



it are small enough to have no effect on the material's 
mechanical behaviour. 

In cellular patterns, especially in dry ones where there 
are no gaps nor overlaps, the movement of each individual 
object is highly correlated with its neighbours'; thus the 
affine assumption is reasonable [T|. In particle assemblies, 
it might apply to dense assemblies of repelling particles, 
which cannot be too close nor too far from each other. 

Whenever this assumption is valid, it considerably sim- 
plifies the description of the pattern evolution. Consider 
for instance a link, £ = r 2 — r\ (eq. [l]). Its time derivative 
is 

— = v (r 2 ) - v (n) . 

Thus, under the affine assumption (eq. [26)) . the velocity 
gradient modifies all links in almost the same way by (see 
eq.l62|): 

d£ affine _ j B ^ 



dt 



Vv { 



Alternatively elasticity and plasticity are defined throuj 
dynamics. A given region of the pattern is said to be in 
elastic, plastic or viscous regime, according to the con- 
tribution to the stress that dominates locally (24]. The 
elastic strain E e i contributes to the reversible part of the 
stress. Plasticity describes the irreversible contribution to 
the stress in the low velocity limit (note that rearranging 
patterns can often deform a lot without breaking). Both 
are solid behaviours, that is, exist in the limit of very low 
velocity gradient. The viscous contribution to stress is ir- 
reversible: it is due to, and thus increases with, the velocity 
gradient; that is, relative movements of objects within the 
material. 



4.3 Link between discrete and continuous descriptions 



In the definition (eq.[7Tbis) of C , the velocity gradient 
h can be taken out of the average: 



That is. 



= Vt>* (£®l). 



affine 

C ~ M Vv. 



(28) 



(29) 



By injecting eq. (j29[) into eq. i|T7|) we show that W is a 
statistical equivalent of the velocity gradient Vv: 



W = M~ X C a ~ e Vv. 



(30) 



This is why we included M 1 only on the left side of W 



4.3.1 Affine assumption 

The affine assumption is analogous to, but much stronger 
than, eq. f22|) . It assumes that the velocity of each indi- 
vidual object is affine too: 

v (r 2 ) a ~" c v (n) + W • (r 2 - ri) . (26) 

An affine flow field and a non-affine flow field are plotted 
on figure fT2l 

In other words, this affine assumption implies that the 
continous velocity gradient has a meaning down to the 
level of individual objects, and that fluctuations around 



4.3.2 Velocity gradient and total strain rate 

By comparing eq. (|23]l with eqs. (|18|30p we identify the 
statistical and dynamical definitions of the total strain 
rate: 

_ „ affine . , . 

V ~ Etot- (31) 

V thus appears as a statistical measurement of the 
symmetrised velocity gradient Etot ■ When only large scale 
measurements are possible, only Etot can be measured. 
However, when the detailed information on links in avail- 
able to perform statistics, measuring V offers several ad- 
vantages. 
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(i) The signal to noise ratio is optimal, in the sense that 
all the local information, and only it, is used. Each link 
acts as a small probe of the local velocity differences: the 
spatial derivative is taken naturally at the places where 
the objects are, not on the larger scale of RVEs. 

(ii) V is intrinsically based on the material's structure. 
It can be defined and measured even if there are only a 
few objects; or if the standard deviation of their veloci- 
ties is large. Averaging over all links provides a statistical 
measurement of the total strain rate. At no point does the 
definition or measurement of V require any affine descrip- 
tion. 

(iii) Physically, we expect V to play a more general 
role than ktot, because it is based on the individual objects 
themselves. For instance, we expect V to be determinant 
in yielding, and thus in the description of plasticity (and 
possibly fl too) [JJ. Similarly, the material's internal dis- 
sipations are probably more closely related to changes in 
the links than to a large scale velocity gradient: this sug- 
gests that the dissipative contribution to the stress arises 
in general from V rather than from ktot- 



4.3.3 Strain, in the elastic regime 

In this section we consider the particular case where the 
material is in the elastic regime. There is no plastic strain 
rate, k p i = 0. Eq. (j25|) becomes simply: 



Ve 



el clastic 



vt 



Etot- 



(32) 



Thus, in the elastic regime, the elastic strain and the 
total strain rate are not independent physical quantities. 
Combining eqs. lf24|) and ([32)1 shows that 



E e l ~ Etot- 



(33) 



More precisely, at least in the linear elastic regime, one 
can identify two quantities: the symmetrised gradient of 
the displacement field, £tot, which is a function of the past 
history of the sample; and the elastic strain e e i, which is a 
function of state. In fact, in elasticity, both quantities are 
considered as equivalent [23] . 

On the other hand, under the affine hypothesis, Ref. 
[TO] for the linear elastic regime (small strains), and Ref. 
[12J for the non-linear elastic regime (large strains), demon- 
strate that: 

affine 

clastic . . 

U ~ etot, (34) 

The demonstration of eq. (|34l) is similar to that for V (eqs. 
[281I3T]) : it uses the same hypotheses, with the additional 
assumption that M and M commute (which is satisfied 
if Mq is isotropic). 

Eqs. I|33|34p show that in the elastic regime: 



strain are possible, £ e i can be measured as e to t- When the 
detailed information on links in available, measuring e e i 
as U offers many advantages, similar to that of V (section 

An acceptable definition of strain must coincide with 
etot in the linear elastic regime. As a consequence, it also 
implies that it is a conjugate of stress: the scalar product 
of stress by an infinitesimal increment of strain equals the 
increment of energy. This is a dynamical constraint on ac- 
ceptable definitions; but it is a weak constraint (especially 
since the conjugate equation is a scalar relation). In itself, 
it is insufficient to define all components of strain. 

There are thus several families of acceptable definitions 
of internal strain [12|I25|; one family contains an infinity 
of acceptable definitions [26]. Some definitions are par- 
ticularly adapted to a discrete pattern's geometry [27] or 
dynamics [28J. 

Here, eq. lfl4|) is a definition of strain which is: (i) one 
of the definitions acceptable in the whole elastic regime, 
even at large strain (eq. [35} where it coincides [12] with 
a true strain [29]; (ii) probably the only definition valid 
outside of the elastic regime [12], when bubbles rearrange 
and move past each other, that is, when the pattern flows: 
the main advantage of eq. l[T4|) is that it does not require 
the detailed knowledge of each object's past displacement. 



4.3.4 Plastic strain rate, in steady flow 

In the more general case, there is a plastic strain rate, 
kpl 7^ 0, and deformations can be strongly non-affme. Eq. 
([32)) does not hold. The current elastic strain e e \ and the 
total strain rate etot are independent physical quantities; 
e e i can no longer be measured as ktot (whether it can be 
measured as U is discussed in section [4.3. 5j) . 

For instance, if the material flows, the displacement 
of an object relatively to its neighbours can be arbitrary 
large; k p i can become much larger than Ve e i/T>t. In the 
extreme cases of steady flows, independent of time, it is 
possible (in absence of advection) that T>e e i/T>t = 0, and 
eq. ([25]) reduces to: 



steady . 
Etot = Epi- 



(36) 



According to eq. lf20|) . a steady flow with a corotational 
derivative that vanishes (meaning no advection nor rota- 
tion effects, see eq. (8lJ), implies that all the geometrical 
strain rate translates into the topological strain rate: 



p steady y. 



Using the identification of eq. (|31j) , we therefore obtain in 
that case: 

affine 
steady 

P * E P i- (37) 



(35) 4.3.5 Complete identification 



Thus the elastic strain e e i can be measured using two dif- The statistical tools U and P are always defined and mea- 
ferent methods. When large scale measurements of total surable, even out of the elastic regime, or out of the steady 
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regime. If we could identify them with e e i and e p i, respec- 
tively, it would make possible to measure the elastic strain 
in all regimes. This is certainly not possible in general, as 
shown by both following counterexamples [30] . 

In granular systems, due to solid friction in the con- 
tacts, irreversible plastic strains appear before the list of 
contacts changes. In solid networks {e.g. solid foams) with 
no topological change, the bond themselves might behave 
plastically, or they might perhaps undergo buckling in- 
stabilities leading to non-reversible stress-strain curves. 
Those are examples in which plasticity occurs before the 
first topological change. 

Conversely, consider a set of rigid cables which re- 
sist tension, but no compression, and tie them together 
at knots to form a redundant, hyperstatic network. Un- 
der given external forces on the knots, some cables will 
be taut, others will dangle and transmit no force. Upon 
changing the forces, the list of taut, tension-carrying ca- 
bles will change. This can be regarded as a topological 
change. The response, which implies displacements and 
strains, is however reversible and might be called elastic. 
Hence a case for which plasticity begins after the first 
topological change. 

This identification might turn possible in some partic- 
ular cases where one can express the stress as a function 
of kinematical quantities. This seems to be the case for 
foams and emulsions OH]- We hope that in these cases, 
statistical measurements can constitute a coherent lan- 
guage to unify the description of elastic, plastic and fluid 
behaviours, as well as facilitate models and tests. 



5 Summary 

In the present paper, we define tools (Table [l]) to extract 
information from a pattern made of discrete objects, sub- 
ject to rearrangements, within a wide class of complex 
materials made of individual constituents such as atoms, 
molecules, bubbles, droplets, cells or solid particles. They 
characterise quantitatively the mutual arrangements of 
these objects, or more precisely the links between neigh- 
bouring objects. 

Their definition, which can flexibly adapt to the ques- 
tions to be answered, is operational. That is, given an 
experimental or simulated pattern, whether in 2D or 3D, 
there is a well defined method to measure them directly as 
statistics on individual constituents (links between neigh- 
bouring sites). This measurement is easy, and requires 
only a few basic operations on a computer: multiplica- 
tion, average, diagonalisation, logarithm. It is robust to 
experimental noise, even if there is a limited number of 
links. 

M, B (or C) and T characterise the current state of 
the pattern, its geometrical changes, and its topological 
rearrangements, respectively. They are explicitly based on 
the pattern's discrete structure. They can be measured 
locally, for instance on a single biological cell, or grain in 
crystals. But they can also be measured as averages over 
a larger region in space, or as time averages. In foams, 
measuring them smoothens out the pattern fluctuations 



due to the discrete nature of bubbles, and evidences the 
underlying behaviour of the foam as a continuous medium. 

Their statistical counterparts U, V (or W) and P, are 
independent of the pattern's discrete length scale. Each of 
them exists and is valid together in elastic, plastic and 
fluid regimes : they unify the description of these three 
mechanical behaviours. They facilitate the comparison be- 
tween experiments, simulations and theories. In at least 
the linear, affine, elastic regime, we suggest how to iden- 
tify them with the quantities which characterise the conti- 
nous mechanics: elastic strain, total strain rate, and plastic 
strain rate, respectively. From a practical point of view, 
this offers the advantage of measuring these continuous 
quantities with an optimal signal to noise ratio, even with 
few discrete objects. On a fundamental side, this provides 
a physical basis to the description of a continuous medium, 
at any local or global scale, by relating it to the individual 
constituents. Moreover, it provides a coherent language 
common to elasticity, plasticity and fluid mechanics. 

The companion paper [l] illustrates most of these points 
on a detailed practical example. 
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A Measurement techniques 

This Appendix, aimed at non-specialists, lists practical 
advices based on our past experience. 



A.l Averaging procedure 

A. 1.1 Weights 

The average of any quantity x is: 

1 ^ 



w x, 



(38) 



where the sum is taken over all links in the averaging 
region. Here w is the weight of the link: for almost all 
links in the averaging region, w = 1; at the boundaries of 
the averaging region, w decreases to zero, different choices 
being possible (section lA.1.2j) . Here we note: 



For instance, the texture is: 



E 1 



M = (m) = 



J2 w 



(39) 



(40) 
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Fig. 13. Examples of averaging procedures: a link's weigth 
vs its position, here the box range is the segment [—a, a]. 
Dots: "all or nothing". Dashes: "proportional". Solid line: 
"coarse grained", here with a hyperbolic tangent profile mir- 
rored around the origin. 

A. 1.2 Choices of weights 

There are at least three main possible choices for the av- 
eraging procedure (Fig. [13]). Once a procedure has been 
selected, it is important to keep consistently the same for 
all measurements. 

The topology is useful for local information, especially 
for a single site: particle or cell. This is the case for in- 
stance when studying the division of a cell [5]. In that 
case, each link is either included or excluded ("all or noth- 
ing", w = or l). One should at least include the links 
between the site of interest and its neighbours (first shell) . 
Statistics are over a few links only, and are easy to com- 
pute, even sometimes by hand. This defines the z-th site's 
texture as a sum over its n,- neighbours, labelled j: 

1 «< 

M i = -Y,m ij . (41) 

"* 3 = 1 

One can also choose to include the second shell (next near- 
est neighbours), third, or even higher. The total pattern's 
texture (eq. [40]) appears as the average of the n S i te site 
textures, weighted by the site's number of links, and with 
a factor 1/2 because each link is counted twice (each link 
belongs to two sites): 

, "site "site 

2 — 1 J — l l — l 

The same definitions (eqs. l4ll and 142]) apply for B or T. 

The geometry is useful for a continuous description, 
typically to measure U{R, t), P(R, t) or V(R, t) in a RVE 
as a function of space and time. This is the case for the 
examples of foam flow which illustrate this paper. The av- 
erage is over all links in a box which respects as much as 
possible the symmetry of the problem: rectangle, annulus. 
The dimensions of the box determine the scale of averag- 
ing. The measurements are performed on boxes at differ- 
ent positions R. The distance between measurements po- 
sitions cannot be more distant than the box size (it would 
leave gaps between boxes), but they can be closer (thus 
boxes overlap) . If the box dimension is much larger than a 
link, one might choose to neglect links which cross the box 
boundary. But in general, links which cross the boundary 
require attention (especially near the box corner) if an au- 
tomatised image analysis is used. An additional choice is 



required. A first possibility ("all or nothing"), computa- 
tionnally simpler, is to look where the link's center lies: if 
the center lies inside the box, the link is assigned to this 
box (w = 1); else, this link is not counted (w = 0) [31]. A 
variant is to assign half of each link (w = 1/2) to the two 
boxes of the two bubble centres it binds [IT]. A second 
possibility ("proportional"), which yields a better preci- 
sion, consists in weighting the link with w equal to the 
fraction of the link which is inside the box (the remaining 
1 — w is outside) [12] . 

The coarse graining is seldom convenient in the prac- 
tical applications considered here. However, theoreticians 
use it [28], especially for the advection term (see Appendix 
IC.ip [15| . A link at position r is counted in a box at po- 
sition R with a weight w(\r — R\). The coarse graining 
function w is a function which is non-increasing, from 
w(0) = 1 to w = 0, and has an integral equal to 1. Its 
width at half height (that is, where w = 1/2) defines the 
scale of coarse graining. It is a continous and differentiable 
function, so that advected links smoothly enter and leave 
the averaging box, without singularity [28] . 



A.2 Choice of M 

Since the reference texture Mq plays almost no physi- 
cal role, its choice is not very important. It depends on 
the problem under consideration, but once its definition 
is chosen, it should be kept consistently. In practice, the 
choice depends on the available information. Here are a 
few possibilities. 

(i) The most favorable case is when Mo can be mea- 
sured directly. In experiment, this is possible when an im- 
age can be chosen as reference, for instance a stress free 
pattern. In simulation (Fig. [3b), this requires to relax the 
stress under prescribed constraints. 

(ii) Mq can be determined theoretically in some cases, 
such as a set of particles which interaction potential is 
known. This occurs in Fig. |(3fr), where the natural refer- 
ence is the honeycomb pattern with a link size: 

e 2 = 2V3A lmk = 2-^=. (43) 

Here Au n k is the area per link, and Mo = ^0-^2/2 ; A = 
3Aii n k is the area per particle or cell (section [A,3.1j) , 

(iii) In cases such as Fig. |TJ), no reference state is 
known in details. If only (£q) is known, we suggest to 
take Mo as isotropic. Although we do not know any fun- 
damental reason for that, it seems to be satisfactory in all 
practical cases we have encountered. From eqs. (0 [HJ it 
writes, in D = 2 or 3 dimensions: 

Mo = ^Id, (44) 

(iv) In some cases, (£0) ls n °t known but we can es- 
timate it. For instance, in a 2D cellular pattern of known 
average area (A) (Fig. [2]), the comparison with hexagons 
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(eq. [43]) suggests to take approximately : 

{%) - , (45) 

and M = ^I 2 /2. 

(v) The most unfavorable case is when even i s 
unknown. A possibility is to take: 

M « A/ D , 

where A is the average ofthe AjS, Af's eigenvalues. Taking 
the arithmetic average, A = J2i K/D, corresponds to the 
assumption that (£ 2 ) is conserved: (^q) « (^ 2 )- Taking 
the geometric average, A = (Ili^j) 1 ^) corresponds to 
the assumption that TrU = 0, which is close to assuming 
that the material is incompressible (see section [3. 1 .4[) . 

A. 3 The case of 2D dry cellular patterns, especially 
foams 

A. 3.1 Number of neighbours 

In 2D dry cellular patterns, the number of neighbours of 
each cell is variable; but its average over the whole pattern 
is always close to 6 neighbours, and thus 6 links, per cell 
|13pi4j. Since each link is shared by two cells, the number 
of links is 3 times the number of cells. This is also true for a 
moderately wet cellular pattern, if neighbours are defined 
on a skeletonized image. It extends to Voronoi/Delaunay 
definition of neighbours for particles. 

Some cells might meet by four ("4-fold vertex"). In 
that case, we recommend to decide that cells which share 
only a vertex should not be considered as neighbours. This 
choice is consistent with the fact that the texture describes 
the cell shape and arrangement (section rA.3,3|) . Moreover, 
this avoids many artefacts when measuring the Tls. 

A. 3. 2 Cell centers versus vertices 

Aubouy et al. [10] chose to describe a cellular pattern 
(such a 2D dry foam) as a network, each site being a vertex 
(that is, a point where three cells meet). Here, we prefer 
to use cell centers, for several reasons. 

(i) First, and most important: the centers move ac- 
cording to the overall velocity field (while vertices have 
a highly fluctating displacement), thus the affine assump- 
tion (eq. [26]) applies. 

(ii) It is more robust, because a cell center is measured 
as an average over several pixels (while a vertex is a single 
pixel, which position might depend on the image analysis 
procedure). 

(iii) This has the advantage of being more general: it 
applies to all other discrete patterns; and even within cel- 
lular patterns, it generalises to wet foams, and to 3D. 

(iv) Finally, the topological rearrangements are well 
characterised (while on the opposite, if £ was defined as 
the vector between two vertices, a Tl occurence would be 
defined as £ a = Id = 0, so that the contribution of a Tl 
to eq. (frT]) would systematically be zero). 



A. 3. 3 Texture and inertia matrices 

Maree et al. [32] propose to measure the shear modulus by 
considering the variation of cell shapes. Each cell's shape is 
characterised by its inertia matrix, (r (g> r) : it looks similar 
to eq. ([2]) , but it is averaged over the position r of pixels 
inside a cell; thus their description is intra- cellular. Ours, 
averaged over the links between the cells, and thus based 
on the shape ofthe overall pattern, is rather inter- cellular. 

In dry cellular patterns, where there are no gaps be- 
tween cells, nor overlaps, the deformation of each cell is 
highly correlated to the global strain; thus, in this case, 
both descriptions coincide and yield approximately the 
same results [5]. 

For instance, note that the shear modulus is the varia- 
tion of elastic stress with respect to infinitesimal variations 
of U. This measurement is robust [3T]ll2j . As mentioned in 
section I3.1.3[ it is not affected if we multiply M and M 
by a same prefactor; and even if we change Mq, see for 
instance eq. lfl4]) . This is why this particular measurement 
gives similar results with both inertia and texture. 

Here, we prefer to use the texture based on cell centers, 
which is more general, for several reasons. 

(i) It also applies to characterise the strain of wet foams 
(where bubbles are round, and thus each bubble's inertia 
is isotropic). 

(ii) It applies to all other discrete patterns, including 
particle assemblies. 

(iii) Centers, rather than shape, are involved in the 
kinematic description, including eqs. I|18|20p . It could in 
principle be possible to define an equivalent of B (and 
even of V) based on inertia matrix, but its physical mean- 
ing is unclear; and it is probably not possible to define an 
equivalent of T (and P) . 

(iv) It extends to more than one cell; while the inertia 
matrix of several cells can be defined, its physical meaning 
is not relevant to the pattern description. 

Note that in the graphical representation of the inertia 
matrix, the ellipse axis lengths are the square root of the 
matrix' eigenvalues [32]. The advantage is that the ellipse 
elongation is the same as that of the actual cell. Here, 
taking the square root of eigenvalues has no physical sig- 
nification for any matrix (except for the texture), so that 
we plot the matrix' eigenvalues themselves (section l2, 2.2ft , 



B Matrices: notations and definitions 

This appendix is aimed at readers who are not familiar 
with the matrices. We list all standard definitions used in 
the text, from the simplest to the most complicated. 

B.l Matrices 

We work here in a space with D = 3 dimensions. A scalar 
is a simple number; a vector is a list of D numbers; a 
matrix is an array of D x D numbers. All these objects 
are tensors, of rank 0, 1 and 2, respectively. In this paper, 
there also appears J, which is a tensor of rank 3 (for which 
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there exists no particular name), that is, an array of D 3 
numbers. 

A matrix A is an array with components Ay, where 
the indices i,j = 1, 2 or 3: 



An A 12 A 13 
A21 A22 A23 
A31 A 32 A33 

Its trace is the sum of its diagonal terms: 



Tr A = A 



11 



A 2 



A 



33 • 



A 4 



(46) 



(47) 

and has the 



Its transposed A 1 has components Ay- — ^ Tl , 
same trace. Any matrix can be rewritten as the sum of its 
symmetric and antisymmetric parts: 

A + A l A -A 1 



A matrix S is said to be symmetric if it is equal to its 
transposed, S = S l , that is, Sy = Sjt (while an anti- 
symmetric matrix is equal to minus its transposed); by 
definition, the symmetric part of A is always symmetric. 
A symmetric matrix can itself be rewritten as an isotropic 
term and a traceless (or deviatoric) term: 



Tr(5) 
D 



I D +Bcv{S), 



(48) 



where Id is the identity matrix in dimension D (appearing 
in eq. [5]). 

Itself, the deviatoric part can be decomposed in diago- 
nal components, called normal differences, and off-diagonal 
ones. 

To summarize, a matrix has in general 9 independent 
components Ay. They can be rewritten as 3 antisym- 
metric ones, namely (A12 — A 2 i)/2, (A 23 — A 32 )/2, and 
(A31 — Ais)/2; and 6 symmetric ones, namely 1 trace An + 
A22 + A33, 2 normal differences An — A22 and A22 — A 33 , 
3 off-diagonal terms (A 12 + A 2 i)/2, (A 23 + A 32 )/2, and 
(A 3 i + Ai 3 )/2. This means that an antisymmetric matrix 
has 3 independent components, a symmetric matrix one 
6, a deviatoric one has 5, an isotropic one has 1. 

The product between matrices is another matrix: 



(49) 



The product between a matrix and a vector is another 
vector: 

(A -a), =J2 A ikak- (50) 
k 

The scalar product between matrices is a number: 



A : B = A lk B kl = Tr (AB* 



(51) 



The ("euclidian") norm of A is a strictly positive number 
defined in any dimension D as: 



1/2 



All = \/A~4 



Y,(AikA kt ) 



(52) 



Many practical applications regard 2D images. Matri- 
cial notations are valid in any dimension D, and it is 
straigthforward to rewrite them in 2D, see also section 
12.2.31 In 2D, a matrix A has in general 4 independent 
components Ay, where i,j = 1 or 2. They can be rewrit- 
ten as 1 antisymmetric one, namely (Ai 2 — A.2i)/2; and 
3 symmetric ones, namely 1 trace An + A22, 1 normal 
difference An — A22, 1 off-diagonal terms (A12 + A2i)/2. 
This means that an antisymmetric matrix has 1 indepen- 
dent component, a symmetric matrix has 3, a deviatoric 
one has 2, an isotropic one has 1. 



B.2 Diagonalisation 

For a symmetric matrix S, there exist three orthogonal 
axes, called S's eigenvectors (from the German word "eigen 11 
meaning "own"), in which S would be diagonal, see for in- 
stance eq. (J4j) - That is, if we used these axes (instead of 
the original ones) to measure the matrix, it would have 
non-zero terms only along its diagonal: 





RSR 1 = diag S = s 2 
s 3 



(53) 



Here R is the matrix of rotation from the original axes to 
the eigenvectors. The three numbers si, S2, s 3 are called 
the matrix' eigenvalues. We label them in order of decreas- 
ing absolute value: | s 1 1 > | S2 1 > | S3 1 . 

They determine many properties of S, including its 
trace and norm: 



Tr(5) 
IISH 



Sl+ S 2 + S3, 



(54) 



If they are non zero, the inverse of S exists, and it is 
diagonal in the same axes as S: 



diag S 



R 




(55) 



If they are strictly positive, the logarithm of S (see eq. 
[El is defined by rotating to the eigenvectors, taking the 
logarithm of the eigenvalue, and rotate back to the original 

axes: 



logS = -R- 



log si 
logs 2 
log s 3 



R. 



(56) 



By construction, logS is symmetric too, and diagonal in 
the same axes as S. 

The literature of mechanics [T91I33] sometimes uses a 
specific definition of shear. It is characterised by a devia- 
toric matrix with two opposite eigenvalues (s± = — s 2 = S) 
and nothing in the third direction (s 3 = 0). Its amplitude 
S is defined as: 



1/2 



S 



V2 



(57) 
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Fig. 14. Ellipses to represent matrices of eigenvalues si and 
S2, with a solid line to represent positive eigenvalues, (a) Two 
positive eigenvalues are represented by a "crossed ellipse". The 
circle represents an isotropic matrix, s 1 = s 2 . (b) When there is 
one positive and one negative eigenvalue, the circle represents 
si = — S2 > 0; the "coffee bean" ellipse is elongated along the 
positive eigenvalue (si > — S2 > 0); the "capsule" ellipse is 
elongated along the negative eigenvalue (— si > S2 > 0). 



In 2D, a matrix can be represented graphically by an 
ellipse, which axes, represented by solid lines, are in the 
directions 9 and 9 + 90°, and have a length s\ and s 2 , 
respectively; the sign of the eigenvalues is labeled specifi- 
cally by plotting a line for a positive eigenvalue, and no line 
for a negative eigenvalue (Fig. [14]). The ellipse anisotropy 
V = (l s i| _ |s2|)/|si/| = 1 - |«2/si| is between 0, for a 
circle, and 1, for an extremely thin ellipse. The ellipse size 
is characterised by | s 1 1 + | S2 1 - 

If S has strictly positive eigenvalues, si > s 2 > 0, it is 
entirely defined by three numbers: first, its trace Tr(£) = 
si + s-2, equal to the ellipse's characteristic size; second, its 
anisotropy rj = l — s 2 /si, equal to that of the ellipse; third, 
the direction 9 of its largest eigenvalue's axis (0° < 9 < 
180°), which is ill-defined for an isotropic matrix (77 close 
to 0). Its determinant s\s 2 is proportional to the ellipse 
area, but it is not used in the present paper. 



B.3 Outer product 

The outer product (or tensor product) of two vectors a, b 
is the matrix of components a$>y. 

(aih a x b 2 aib 3 \ 
a 2 b\ a 2 b 2 a 2 b 3 . (58) 
a 3 bi a 3 b 2 a 3 b 3 J 

Its trace is the scalar product a ■ b: 

Tr(a <g> b) = a ■ b = a\b\ + 0262 + a 3 b 3 . 

For instance, eq. (JSJ, and eq. ([3} or ((6|, write: 

m = Z®Z, 

M = (£ ® i) . (59) 



In eq. I|59|) , the physical interpretation of the outer prod- 
uct by i is that it transforms the surface integral of a 
discrete individual vector, the link, into a bulk integral of 
a continuous average matrix, the texture [!0lfT2] . 



The outer product is also used for the notation V = 
(d/dri) = (d/dx,d/dy,d/dz) (called gradient, or nabla) 
which symbolises the space derivatives: 



Vv = V ' 



dvi dv2 dv3 

dx ffx ffx 

yvi 0V2 0V3 

dy dy dy 

dv\ dv2 dvs 

dz dz dz 



(60) 



More precisely, eq. 11221) uses its transposed: 



On 



(61) 



Similarly, the rotational is the vector product V» x vj. 

However, rheologists [191133] often prefers the notation 
gradu = dvi/drj. This creates an ambiguity with eqs. 
(|60|6ip . In case of doubt, it is safe to come back to indices, 
which are unambiguous. For instance, the demonstration 
of eq. I|27p writes as follows: 



j 3 



(62) 



In section 12.3.11 and in eq. l(66|) . the notation J = 
v Cg) M means Jijk = ViMjk, and the notation V ■ J 
is equivalent to J^i dJijk/dn. 



C Time evolution 

C.l Finite time interval At 

This appendix can be useful for a user who wants to anal- 
yse a movie (and not only a static image). It helps to 
understand the definitions, units and measurements of B 
and T; as well as the time evolution of M (eq. [9| . 



C.l.l Time interval between images 

We consider a movie. To obtain good measurements, av- 
erages can be performed on a large time interval r, that 
is, a large number r j At of images. Here At is the time 
interval between two consecutive images, at times t and 
t+At. Ideally, measurements should not depend too much 
on the exact value of At. In practice, however small At 
is, it is finite, and this discretisation has consequences, see 
below eqs. I)64|68j) 

If possible, At should be chosen small enough to en- 
able a good tracking of objects from one image to the 
next (typically, during At, relative displacements of sites 
should be a fraction of the inter-sites distance). This is 
unimportant if objects are labelled individually, as is the 
case in simulations. 

At should also be chosen small enough that for any 
time-dependent (but space- independent) quantity x, Ax/ At 
tends towards its time derivative dx/dt. When x depends 
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both on space and time, Ax /At tends towards dx/dt if 
the measurement box moves along with the links (so-called 
"Lagrangian" point of view, often useful in theory). When 
the measurements are performed in a fixed region of space 
(so-called "Eulerian" point of view, often useful in prac- 
tice, especially in steady flows), Ax/ At tends towards its 
partial time derivative dx/dt. 
Eq. (|40|) can be rewritten as: 



N tot M = J2 w£ ® £ = J2'' 



(63) 



Its variation between successive images, divided by At, in- 
volves the links t a (resp. td) appeared (resp. disappeared) 
during At. This means that N to t and M are evaluated 
both at t and t + At; t a and m a (resp. id and ma) are 
evaluated at t + At (resp. t). Then: 



Am 



A{t®t) 
A£®£ + £< 



At + £At <g> At. 



(64) 



Here £ = 1 if t is evaluated at time t; £ = — 1 if £ is 
evaluated at time t + At; and £ = if £ is the average 
of its values at t and t + At, which is recommended to 
simplify eq. (|68|) . 



C.1.2 Effect of discrete time 



The variation of eq. ([63|l between successive images writes: 



AjNtotM) 
At 



EAw 
— m 



E 



Am 



w- 



At ^ At 

J2 a w a m a - J2d w dm d 
At 



(65) 



Quantities in eq. I|65| are extensive, and are convenient 
for practical measurements. However, for theory, intensive 
quantities are easier to manipulate (Table [1]). Dividing 
both sides by Ntot yields: 



AM M AN t , 



At 



N. 



At 



-V -Jm + B + T. 



(66) 



Eq. (|66j) tends towards eq. ([9]) in the limit of small At, as 
we now show term by term. 

The relative variation of Ntot during At is A log Ntot/ At, 
and at small At it tends towards d log N to t/dt. 

The advection term —V • Jm is the term in Aw/ At, 
in eq. f65|) . It is due to links entering or exiting the region 
where M is measured (see section [2.3. ip . In a first ap- 
proximation, Jm — v®M. The demonstration is delicate 
and we do not develop it here. Briefly, when the averag- 
ing procedure uses a coarse-graining function w{r{t)), it 
is possible to transform a time derivative of w into a space 
derivative; this involves dr/dt, that is, the local velocity 
[151128] . 

The geometrical variation term is: 



B = 



N r 



N,, 



Am 
~At 



(67) 



Here iV c is the number of links conserved between both 
images, that is, the number of terms involved in the aver- 
age noted (.). Most links contribute to eq. I|66p . but each 
one has a small contribution. If At is small enough, then 
At I At tends towards dt/dt and the correction N c /N to t 
tends towards 1, so that we obtain eq. ljlO| . Eq. (|64| . im- 
plies that: 

B = C + C l +O(0, (68) 
where C is defined as: 



C = 



Nr. 



N t , 



!<g) 



At 
At 



(69) 



If we choose £ = 0, eq. (JHHJ) means that B is twice the 
symmetrical part of C: 



B = C + C t + 0{i), 
If At is small enough, eq. (|69|) yields: 



(70) 



l dt 



■{xm (vm izmy 

<*#> <r#> Mf) 



(71) 



In general, C is not symmetric; its transposed is C = 
(dt/dt®t). 

The topological term is between the last parentheses 
in eq. (|65| . Note that we can derive its exact prefactor: 



r = — <m a ) mj 



(72) 



where ANd is the number of disappeared links (that exist 
at t but no longer at t + Z\i); and ziiV a the number of 
appeared links (that exist at t + At but not yet at t) . In 
the limit of small At, eq. lf72j) tends towards eq. (fTTj) . In 
eq. (|72p . the number of terms involved in the average (m a ) 
(resp: (m^)) is /A^V a (resp: ANd), which is much smaller 
than Ntot- This represents a small number of links, each 
having a large contribution to eq. (|66|l . Thus statistics on 
T are always much noisier than that on M or B, and it 
is advisable to integrate T over a long time r. 



C.2 Objective time derivatives 

This technical section is rather aimed at specialists. It dis- 
cusses the different objective derivatives which appear in 
the course of this paper, when estimating the time deriva- 
tives of M or U. By definition, the objective derivative of 
a matrix is invariant after a change in any other rotating 
frame of reference. It expresses the fact that M or U are 
intrinsic properties of the material. In practice the objec- 
tive derivative writes as the sum of a total (Lagrangian) 
derivative plus corrections involving the velocity gradient. 
In principle there is an infinity of possible objective deriva- 
tives. The following sections show that objective deriva- 
tives can be selected and calculated for the evolution of 
matrices M and U. 



20 



F. Graner et al.: Discrete rearranging patterns: Robust tools 



C.2.1 Time evolution of M 



C.2.3 Kinematic equation of U 



In the affine assumption, the geometrical term B of eq. 
^ can be rewritten using eqs. (|17|30|70[) : 



B = MW + WtM 



affinc 



MVv + Vv l M. 



(73) 



Thus, the second term of the right hand side of eq. |(9|) 
can be grouped with its left hand side, formally appear- 
ing as a Maxwell upper convective tensor derivative (see 
e.g. [34]): 



v 
M 



dM 



VM-Vd* M MVv. 



(74) 



Here we have assumed incompressibility to transform V • 
J M w V • (v g) M) into v ■ VM. Hence eq. © appears 
as a conservation equation for M: 



In section lC.2.2( with U small and M isotropic, M com- 
muted with its time derivative. Thus, knowing the time 
derivative of M, eq. (f!4|) immediately yielded the time 
derivative of U. This enabled to eliminate M from the 
time evolution of U (eq. [78]) . 

This is not the case in general. If the eigenvectors of M 
change (rotate) with time, M does not necessarily com- 
mute with its time derivative. There is no simple relation 
between the time derivatives of M and U. 

It is thus tedious to obtain the time evolution of U 
versus U, instead of versus M. We do not develop here the 
calculation. Briefly, we start from eq. |77j) . We differentiate 
it: 



dM 



Mn 



' dt 



O [U 2 



(79) 



We then inject eq. |79j) in the time evolution of M (eq. 
[9]) at lowest order terms in U ; then eliminate M using eq. 
(El: 



M 



all in 



(75) 



which source T is due to topological changes. Note that 
the present approach has unambiguously selected the up- 
per (rather than the lower, or any other) convective tensor 
derivative. 



C.2.2 The small U assumption 



Inverting eq. (fl"4]) . the texture develops as: 
M = M cxp (217) . 



(76) 



In a plastic material such as considered here, the elastic 
internal strain is seldom much larger than unity. This is 
the case for foams, where deformation of bubbles does not 
excess the size of two bubbles before topological changes 
occurs, and thus \\U\\ is bounded. In the case where the 
strain U is small everywhere, eq. f76j) becomes: 



M = M [I + 2U + O 



[U 2 )] s ~ n M 



2M U. (77) 



dU 

^- = -v-ju + v-nu 

dt 

Ufl 1 +0 (f2U 2 ) - P. 



(80) 



The term —P appears on the r.h.s. of eq. (|80| thanks 
to the factor —1/2 in eq. (J20j) - Physically, we can track it 
across section [C.2.2[ back to the factor 2 in eq. I|76p . thus 
in eq. (fl4| . 

The higher order terms 0(f2U 2 ) are often negligible 
in a plastic material such as considered here, where the 
elastic internal strain is seldom much larger than unity. 

On the other hand, the advection term and the 
rotation term flU + Ufl 1 have symmetries which are dif- 
ferent from that of V and P. They may thus not neces- 
sarily be negligible. They can be regrouped using the total 
corotationnal ("Jaumann" [35]) objective derivative: 



VU 
~Vt 



dU 
Iff 



— = ^- + V-Ju + f2U- Uf2, 



(81) 



where we recall that fl l — —17. We thus approximately 
obtain eq. (|2"T]) . This result is close to the objective deriva- 
tive for the logarithmic Hencky strain (up to a small cor- 
rection on the rotation rate), see 1361137] . 



Let us also assume that the reference configuration is 
isotropic, Mq = MqI. Then, eq. (|75| becomes: 



2M I U - t 



small 



or equivalently: 



u s 1 V P, 



(78) 



where P = — T/(2M ) is the rate of the plastic strain. 
This linear assumption thus simplifies eqs. (|21|8ip . Again, 
this selects an upper convective tensor derivative. 
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